Matrix Operations Demonstration

With the iPython notebook and pylab support, we have a very nice interactive shell for scientific computing.

The ipython notebook lets you mix rich text including equations, computation, and graphics into a single working document. You can easily revise computational steps without having to re-run an entire script. Notebooks can be saved and reloaded, printed, and accessed remotely.

Example equation: $X = \sum_{i=1}^{n}i^2$ produces: $X = \sum_{i=1}^{n}i^2$

Pylab is the name given to a common collection of python scientific computation and plotting libraries, imported into a common namespace to be more user friendly when used interactively. Pylab imports numpy and matplotlib, a MATLAB-like plotting library.

Invoked via the command line as: ipython notebook --pylab inline

Rich text is provided using the "markdown" syntax. For information on markdown, see http://en.wikipedia.org/wiki/Markdown.



In [1]:

    
A = zeros((4,4))
I = eye(4)
print "A=",A, "\nI=", I









    



A= [[ 0.  0.  0.  0.]
 [ 0.  0.  0.  0.]
 [ 0.  0.  0.  0.]
 [ 0.  0.  0.  0.]] 
I= [[ 1.  0.  0.  0.]
 [ 0.  1.  0.  0.]
 [ 0.  0.  1.  0.]
 [ 0.  0.  0.  1.]]

When using pylab, many common MATLAB-like functions exist in our default namespace. A good resource for those coming to python from MATLAB is http://www.scipy.org/NumPy_for_Matlab_Users.



In [2]:

    
N = randn(4,4)
N









    Out[2]:





array([[-1.09982708, -0.66931755, -0.4632849 , -0.9946823 ],
       [-0.27885558,  1.44037606, -0.0579929 , -0.26660722],
       [ 0.38743295,  1.18197318,  0.70542364,  0.13824228],
       [ 0.05602669, -0.53328956,  1.1624428 , -0.22843116]])

Note the data type is a numpy multidimensional array. </br> Addition and scalar multiplication on arrays work as expected.



In [3]:

    
N+4*I









    Out[3]:





array([[ 2.90017292, -0.66931755, -0.4632849 , -0.9946823 ],
       [-0.27885558,  5.44037606, -0.0579929 , -0.26660722],
       [ 0.38743295,  1.18197318,  4.70542364,  0.13824228],
       [ 0.05602669, -0.53328956,  1.1624428 ,  3.77156884]])

The * operator between arrays performs component-wise multiplication. For matrix multiplication, use dot(X,Y) instead.



In [4]:

    
print N*I #array component mult.
print 
print dot(N,I) #matrix mult. using arrays









    



[[-1.09982708 -0.         -0.         -0.        ]
 [-0.          1.44037606 -0.         -0.        ]
 [ 0.          0.          0.70542364  0.        ]
 [ 0.         -0.          0.         -0.22843116]]

[[-1.09982708 -0.66931755 -0.4632849  -0.9946823 ]
 [-0.27885558  1.44037606 -0.0579929  -0.26660722]
 [ 0.38743295  1.18197318  0.70542364  0.13824228]
 [ 0.05602669 -0.53328956  1.1624428  -0.22843116]]

There is also a matrix data type. </br> The multiplication operator * performs matrix multiplication, and the multiply(X,Y) performs component-wise multiplication.



In [5]:

    
N2 = matrix(N)
I2 = matrix(I)
print N2 * I2 #matrix mult
print
print multiply(N2,I2) #component mult









    



[[-1.09982708 -0.66931755 -0.4632849  -0.9946823 ]
 [-0.27885558  1.44037606 -0.0579929  -0.26660722]
 [ 0.38743295  1.18197318  0.70542364  0.13824228]
 [ 0.05602669 -0.53328956  1.1624428  -0.22843116]]

[[-1.09982708 -0.         -0.         -0.        ]
 [-0.          1.44037606 -0.         -0.        ]
 [ 0.          0.          0.70542364  0.        ]
 [ 0.         -0.          0.         -0.22843116]]

With PyLab, we can easily call the scipy linear algebra functions normally found in scipy.linalg package.</br> (See: http://docs.scipy.org/doc/scipy/reference/linalg.html)



In [6]:

    
(evs, evecs) = eig(N) #if input is ndarray, output is too.
(evs2, evecs2) = eig(N2)  #if input is matrix, output is matrix
print type(evecs), type(evecs2)









    



<type 'numpy.ndarray'> <class 'numpy.matrixlib.defmatrix.matrix'>



In [7]:

    
#if not using pylab, do the following
import scipy.linalg as LA 
(evs, evecs) = LA.eig(N)

We can do some easy MATLAB-like plotting too!



In [8]:

    
X = randn(1000)
hist(X)









    Out[8]:





(array([  5,  25,  51, 138, 235, 230, 180, 100,  28,   8]),
 array([-3.15868084, -2.5495798 , -1.94047876, -1.33137772, -0.72227668,
       -0.11317564,  0.4959254 ,  1.10502644,  1.71412748,  2.32322852,
        2.93232956]),
 <a list of 10 Patch objects>)



In [9]:

    
Y = 3*X
plot(X,Y)









    Out[9]:





[<matplotlib.lines.Line2D at 0x10c1d7d50>]



In [9]: